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1 Introduction 



Jet physics plays an important role at the TEVATRON and will become even more important 
at the LHC. It provides information on the strong interactions and forms quite often important 
backgrounds for searches of new physics. While jet observables can rather easily be modelled at 
leading order (LO) in perturbation theory [1-10], this description suffers from several drawbacks. 
A leading order calculation depends strongly on the renormalization scale and can therefore 
give only an order-of-magnitude-estimate on absolute rates. Secondly, at leading order a jet 
is modelled by a single parton. This is a very crude approximation and oversimplifies inter- 
and intra-jet correlations. The situation is improved by including higher-order corrections in 
perturbation theory. 

At present, there are many NLO calculation for 2 — > 2 processes at hadron colliders, but 
only a few for 2^3 processes. Fully differential numerical programs exist for example for 
pp^3jets [11-13], pp + 2jets [14], pp -> ttH [15,16] and pp^// + 2jets [17,18]. 
The NLO calculation for pp — > tt+jet is in progress [19]. In the examples cited above the 
relevant one-loop amplitudes were usually calculated in an hand-crafted way by a mixture of 
analytical and numerical methods. However it has become clear, that this traditional way reaches 
its limits when the number of external particles increases. On the other hand, it is desirable 
to have NLO calculations for 2 — > n processes in hadron-hadron collisions with n in the range 
of n = 3, 4, 6, 7. QCD processes like pp — > n jets form often important backgrounds for the 
searches of signals of new physics. To overcome the computational limitation, there were in the 
last years several proposals for the automated computation of one-loop amplitudes [20-33]. 

In this paper we report on the implementation of an algorithm for the automated computation 
of one-loop integrals, which occur in a massless quantum field theory into a numerical program. 
For QCD processes at high-energy colliders the massless approximation is justified for all quarks 
except the top quark. The number of external particles of the loop integrals is not restricted within 
our approach. All calculations are done within dimensional regularization. When combined with 
the appropriate contributions coming from the emission of an additional parton, the project we 
report on here will provide a numerical program for the automated computation of 2 — > n NLO 
processes in massless QCD. As our approach is valid (in theory) for all n, the actual limitation 
on n will result from the available computer power for the Monte Carlo integration. 

The problem which we address in this paper is the fast and efficient numerical evaluation of 
scalar and tensor one-loop integrals in a massless quantum field theory. Tensor integrals are loop 
integrals, where the loop momentum also appears in the numerator. Loop integrals are classified 
according to the number n of internal propagators (or equivalently the number of external legs), 
as well as the rank r, counting the power to which the loop momentum occurs in the numerator. 
It is a well known fact, that all one-loop integrals can be expressed in terms of the scalar two-, 
three- and four-point functions, up to some trivial extra integrals, which are mainly related to a 
specific choice of the regularization scheme. The task is to calculate numerically the coefficients 
in front of the basic scalar integrals. It is tempting to do this with a single algorithm, which covers 
all cases in a uniform way. Although several of these algorithms exist, a particular algorithm will 
perform well for most configurations, but can lead to numerical instabilities in certain corners of 
configuration space. We therefore opted for a "patch-work"-style, treating loop integrals with n 
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propagators and rank r on an individual basis. This reduces to a certain extent the dependency 
on the caveats of a particular algorithm and allows us rather easily to replace in future releases 
of the program a particular reduction method with an improved version. 

We employed the following strategies for the reduction of one-loop integrals: The two-point 
functions are rather easy and are therefore evaluated directly. For the reduction of tensor integrals 
with n > 3 we use spinor methods and follow mainly the recent work by del Aguila and Pittau 
[29, 34-36]. This leads to scalar integrals, where additional powers of the £-components of the 
loop momentum can still be present in the numerator. If such powers are present, the resulting 
integrals are rather easy and are evaluated directly. It remains to treat scalar n-point integrals with 
n > 5 and to reduce them to the basic set. For n = 5 and n = 6 the reduction is unique [25,37-39]. 
This is no longer true for n > 7. In the latter case we use a method based on the singular value 
decomposition of the Gram matrix [27,40]. These steps reduce all integrals to the basic set of 
scalar two-, three- and four-point functions. The latter are then evaluated in terms of logarithms 
and dilogarithms. 

This paper is organized as follows: In the next section we introduce our notation. Section 
3 discusses the reduction of tensor integrals. Section 4 evaluates higher-dimensional integrals, 
resulting from additional powers of the £-components of the loop momentum in the numerator. 
Section 5 treats the reduction of scalar n -point integrals for n > 5. In section 6 we comment on 
the numerical implementation. Finally, section 7 contains our conclusions. In an appendix we 
provide the necessary details on spinors as well as the explicit expressions for the basic scalar 
integrals and methods for the numerical evaluation of some special functions. 

2 Definitions and conventions 

The general convention for a scalar one-loop n-point integral is 
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d D k 
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(1) 



m§ k 2 (k 



Pl) 2 ...(k-pi - ...pn-l) 2 ' 



with D = 4 — 2e. We further use the notation 



H=Y<Pji ki = k-qt. 

7=1 

The flow of momentum is shown in fig. (1). The kinematical matrix S is defined by 



(2) 



Sij = (qi-qj) , 



(3) 



and the Gram matrix is defined by 



Gij = Iqiqj. 



(4) 



Integrals of the type 



I 



(5) 



3 



Pi Pn 

Figure 1 : The labelling for a generic one-loop integral. The arrows denote the momentum flow. 



are called tensor integrals. These integrals are said to have rank r, if the loop momentum appears 
r-times in the numerator. These integrals are always contracted with a coefficient J£ v „ r , which is 
a product of n tree-level currents. This coefficient depends on the momenta and the polarization 
vectors of the external particles of the scattering process. Since trees can be attached to the 
external lines of a one-loop integral, the external momenta pj of a one-loop integral are in general 
not the momenta of the external particles in the scattering process, but rather sums of the latter. 
The coefficient 7" ^ can be computed efficiently in four dimensions. 

It is therefore appropriate to discuss different variants of dimensional regularization. The 
most commonly used schemes are the conventional dimensional regularization scheme (CDR) 
[41], where all momenta and all polarization vectors are taken to be in D dimensions, the 't 
Hooft-Veltman scheme (HV) [42], where the momenta and the helicities of the unobserved par- 
ticles are D dimensional, whereas the momenta and the helicities of the observed particles are 
4 dimensional, and the four-dimensional helicity scheme (FD) [43-45], where all polarization 
vectors are kept in four dimensions, as well as the momenta of the observed particles. Only the 
momenta of the unobserved particles are continued to D dimensions. 

The conventional scheme is mostly used for an analytical calculation of the interference of 
a one-loop amplitude with the Born amplitude by using polarization sums corresponding to D 
dimensions. For the calculation of one-loop helicity amplitudes the 't Hooft-Veltman scheme and 
the four-dimensional helicity scheme are possible choices. All schemes have in common, that the 
propagators appearing in the denominator of the loop-integrals are continued to D dimensions. 
They differ how they treat the algebraic part in the numerator. In the 't Hooft-Veltman scheme the 
algebraic part is treated in D dimensions, whereas in the FD scheme the algebraic part is treated 
in four dimensions. It is possible to relate results obtained in one scheme to another scheme, 
using simple and universal transition formulae [46-48]. 

Since the efficient numerical calculation of the coefficient J" „ relies on the Fierz identity 
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in four dimensions, we are lead to the choice of the four-dimensional helicity scheme. In this 
scheme we can assume without loss of generality that the coefficient 7^ is given by 

•C.^ = <«i-|Y w |^i->-<«r-|Y Pr |^->, (6) 

where (a; — | and \bj— ) are Weyl spinors of definite helicity. It is convenient to denote spinor 
inner products as follows: 

(pq) = (p-\q+), [qp] = (q+\p-)- O) 
Important relations satisfied by the Weyl spinors are 

(p-Mq-) = (q+Mp+), (8) 

and the Fierz identity 

(a-\y II \b-)(c+\f\d+) = 2(ad)[cb]. (9) 
Therefore we consider tensor integrals of the form 

f = e^u 2 Ha,-\y \bi-) (a-\y \b-)[— ^'"^ 

(10) 

where hf^ denotes the projection of the D dimensional vector W onto the four-dimensional sub- 
space. A peculiarity of the four-dimensional helicity scheme is given by the fact that the dot 
product of kfL, with itself does not cancel exactly a propagator, i.e. the D-dimensional k 2 D ^ is 

given as the sum of the four-dimensional k 2 ^ and k 2 _ 2z y consisting of the £-components: 

k\ D) = k 2 4) +k 2 _ 2e) . (11) 

When no conflicting interpretations are possible, we will often drop the indication of the dimen- 
sion of the underlying space. As a consequence we have to consider a generalization of eq. (10) 
by allowing additional powers of ^_ 2e ) m me numerator: 

7" - -\y \bi-) (a-\y \b-)[^ H-2 £ )) ^)-% 

(12) 

The result of eq. (12) can be expressed in the form 

I"/ = ^+%i+Co+0(8). (13) 
£ z £ 

We are mainly interested in the coefficient Co. Besides that, the knowledge of the coefficients 
C 2 and C_i provides additional cross checks, as the divergent part of the Laurent series has to 
cancel against similar parts coming from the real emission and renormalization. The purpose of 
this paper is to set up a scheme for the numerical calculation of the coefficients C_2, C_i and Cq. 
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3 Tensor reduction 



The classical method for the reduction of tensor one-loop integrals is the Passarino-Veltman 
algorithm [49-52]. Here, we use instead spinor methods, discussed for example in [29,34-36]. 
The spinor methods have the advantage that they avoid to a large extent the appearance of Gram 
determinants, or in cases where they cannot be avoided, reduce them to square roots of Gram 
determinants. Alternative methods, like for example approaches based on dual vectors or raising 
and lowering operators are discussed in [53-57]. In this section we give an algorithm for the 
reduction of integrals of the form as in eq. (12) towards integrals of the form 

7°* = e^u 2 *[— ("*(-*)) (14) 

" ^ J iiS k 2 (k-pi) 2 ...(k-pi-... Pn -i) 2 ' U4j 

We do this by treating the different cases of n separately: Two-point functions (n = 2) are rather 
simple and are calculated directly in section 3.2. For n > 3 we use spinor methods. The cases 
n = 3 and n = 4 are special, as there are only two, respectively three independent external mo- 
menta. The tensor three-point functions are discussed in section 3.2, while the tensor four-point 
functions are treated in section 3.4. Finally, for n > 5 we use a general method for the tensor 
reduction, which is discussed in section 3.5. 



3.1 Generalities 

The basic formula for the Passarino-Veltman algorithm states, that the scalar product of a loop 
momentum with an external momentum reduces the rank of the tensor integral: 

2 Pi -kj = kf^-kf + q 2 -qf^-lpiqj. (15) 

This formula is valid independently of which variant of dimensional regularization is used, since 
the £-components cancel between k 2 _ x and k 2 . Therefore the subscript indicating if the loop 
momentum kj lives in D or four dimensions was dropped. 

The first step for the construction of the reduction algorithm based on spinor methods is to 
associate to each n-point loop integral a pair of two light-like momenta l\ and h, which are linear 
combinations of two external momenta p, and pj of the loop integral under consideration [29]. 
Note that pi and pj need not be light-like. Obviously, this construction only makes sense for 
three-point integrals and beyond, as for two-point integrals there is only one independent external 
momentum. If p, and pj are light-like, the construction of l\ and h is trivial: 

h=Pi, h = Pj- (16) 
If pi is light-like, but pj is massive one has 

h=Pu h = -U2Pi + Pj, (IV) 



where 



cc 2 = (18) 
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The inverse formula is given by 

Pi = h, Pj = a 2 h+h- (19) 
If both pi and pj are massive, one has 

/l = l — ]rz-(Pi- a iPj)> h = z — X — -(-a 2 pi + Pj). (20) 
1 — aia2 1 — aia2 

If Ipipj > 0, ai and a 2 are given by 



2pip:-VA 2p i p j -y r K 

- J —? i , 0C2 = —r 



ai = FU i 2 , a 2 = , J 7 ■ (21) 



If IpiPj < 0, we have the formulae 

_ 2p iP j + y/A _ 2p lPj + VA 

«1 r-j , «2 r-2 • (22) 

2/?j 2/?f 

Here, 

A = (Ipipjf-Arfp). (23) 

The signs are chosen in such away that the light-like limit pj^O (or pj — > 0) is approached 
smoothly. The inverse formula is given by 

p i = k + a l l 2 , pj = a 2 h+l 2 . (24) 

Note that l\, l 2 are real for A > 0. For A < 0, l\ and l 2 acquire imaginary parts. These formulae 
can be used in the following ways: First we may decompose any four- vector p into a sum of two 
null- vectors: 

p = an + l, (25) 

where n is an arbitrary null-vector and 

P 2 

I = -an + p, a=-^—. (26) 
2pn 

Secondly, we may decompose ^ as follows: 

1 



2/1/2 



[(2kil 2 ) h + (2kih)h -hHih ~ HHih\ , (27) 



where l\ and l 2 are obtained from decomposing p, and pj into null- vectors. Note that this for- 
mula can be proved by solely using the anti-commutation relations for the Dirac matrices and is 
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therefore valid in the HV/CDR-scheme as well as in the FD-scheme. The main application for 
eq. (27) will be the application towards the spinor strings 

{a-\li {A) \b-) = ^[(2kh)(a-\h\b-) + {2kh){a-\h\b-) 

- (ah) [hb] (h ~ \H{A) | /i-> " (ah) [lib] (h ~ \H{A) \ h~)] , (28) 

appearing in eq. (12). We note that the scalar products of ki with h or h are linear combinations 
of the scalar products of k[ with pi and pj, 

2kh = - — [2pik - ailpjh] , Ikih = - — [-azlpiki + 2 Pj ki] , (29) 

1 — a\U2 1 — aia2 

and therefore immediately reduce the rank of the tensor integral through eq. (15). Eq. (28) 
allows us to replace an arbitrary sandwich 

(a-\y [4) \b-) (30) 

with the standard types 

</l-|^ ( 4)|/2-> and (h-\H(4) \h~), (3D 

plus additional reduced integrals. This procedure can easily be iterated. Note that l\ and h 
depend on the external momenta of the loop integral. In general, these two vectors have to be 
re-defined when pinching a propagator. 

3.2 The tensor two-point function 

The tensor two-point function is special, as it does not fit into the general scheme, which we use 
for the tensor reduction. This is due to the fact that the two-point function depends only on one 
external momentum. Fortunately, the two-point function is simple enough, such that one can 
solve the problem by direct calculation. We consider the general tensor two-point integral 

h ~ e ^ J fat I (k-p)W ( } 

1 D 

= e £ ^p 2£ Jda J ^ (^-kf_ 2£) y (k + apf l ...(k + apf r [-k 2 + a(l -a) (-p 2 )Y\ 
Expanding (k + ap) m ... (k + ap) Mr yields terms of the form 

a r-2tj^ a(l) _ _ j^% (20 pMc(2 t +l) _ _ p^(r) _ (33) 

Note that terms with an odd number of k^'s vanish after integration. We further have 

d -4^ f (k 2 ) = g — d -4k 2 f(k\ 

J il& u J in? 

^Ak^Wk^Hk 2 ) = 1 ^ + \ g ^ ^(k 2 ) 2 f(k 2 ). 

J m§ D(D + 2) J K ' JK ' 
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In general we have 

r d D k r(-) ( d D k 

rL^k^...k^f(k 2 ) = 2- w ,,} 2 \ (g^ W2 ■ ■ .g^-w™ + permutations) / (k 2 ) w f(k 2 ). 

J i%? Y \^+ w ) J tod 

The fully symmetric tensor structure 

5 w--Ww = gtnK g Ww-iKw + permutations (34) 

has (2w — 1)!! = (2w — l)(2w — 3)...l terms. We obtain in the absence of powers of k 2 _ 2e ) 



^ f daa r - 2t f ^-k»K..kK* [-k 2 + a(\-a) {-p 2 )}~ 2 = 
J J ill 2 



o 

P 2 \ 1 s ^ 2t r(l + r-f-e)r(2-2e) ^ 



2 y r(l-e)r(2 + r-2e) 

= (-f ) ^-■■■^|^{l+ e [2Z 1 (r+l)-Z 1 (r-0-2] + o(8 2 )}7 2 , (35) 

where Zi (n) is a harmonic sum 



Zi(n) = (36) 

7=1 J 



and 7 2 is the scalar two-point function: 

Since h starts at 1/e we can neglect (e 2 ) terms in eq. (35). If powers of £ 2 _ 2e ) are present, we 
obtain if all indices are contracted into four-dimensional quantities 

l D 

e e V e / daa r - 2t J " (-fc 2 _ 2e) )V..Jfc«' [-k 2 + a(\ -a) {-p 2 )]~ 2 = 



.2 



t 



(s-\)\{r + s-t)\ 



- -e(P 2 )^-%j^--- ^ ( -;; )! -7 2+0( e) 

+terms, which vanish when contracted into 4-dimensional quantities. (38) 



3.3 The tensor three-point function 

For tensor three-point integrals we may use eq. (28). The first two terms on the r.h.s. of eq. (28) 
reduce the rank immediately. We can therefore assume that the tensor structure is a product of 

(l x -\^\l 2 -) and (Z 2 -|#j%-). (39) 
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Note that the index of the loop momentum is irrelevant, 



(h - W^lh-} = (h - W^lh-} = (h - 14 4) |/ 2 ->, (40) 

since the following sandwiches vanish: 

(h-Wi\h-) = (h-W2\h-)=0. (41) 
Two different spinor types reduce the rank: 

(/i-|^ (4) |/2-)(/ 2 -|^ (4) |/i-) = (2/i^)(2/ 2 fc z )-(2/ 1 / 2 )(fc / (4) ) . (42) 

Here, (kf^ denotes the square of the four-dimensional components and does not exactly cancel 
a propagator. Using 

*(4) = k2 (D)-k 2 { - 2£ ) (43) 

together with eq. (56) will lead to integrals in higher dimensions. It remains to discuss the case 
of a tensor structure of the same spinor type, e.g. either 

(/ 1 -|^ (4) |/ 2 -)...(/ 1 -|^ 4) |/ 2 -), (44) 

or the same situation with l\ and / 2 exchanged. It is easy to see that these terms will vanish after 
integration, since any contraction of 

</i-|Y w |fa->...(ii-|Yivlfa-> (45) 

with pj, p 1 ^ or g^ v will vanish. 



3.4 The tensor four-point function 

For the four-point function two new features appear: One can no longer shift freely the loop 
momentum inside the spinor sandwiches and tensor structures of the same spinor type, as in 
eq. (44), no longer vanish identically. On the other hand, the four-point function has, apart 
from the two external momenta p, and pj used to construct l\ and / 2 , one additional independent 
external momentum, labelled in the following. For the tensor reduction one starts again with 
eq. (28), possibly preceded by a shift in the loop momentum, which synchronizes all occuring 
loop momenta in the numerator from 

to *f ...*<*. (46) 
It is therefore sufficient to consider a tensor structure, which is a product of 

(/i-l^lfe-) and {l 2 -\^\l x -). (47) 
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If in the tensor structure both spinor types appear, we can use eq. (42): 

(h - WPlh-Hh-WPlh-) = (2Zi*0 (2/2*/) - (2/i/ 2 ) (*{ 4) ) 2 - 

If on the other hand in the tensor structure only a single spinor type occurs, we now use the third 
external momentum p 3 and write: 

(/1 - i^i/ 2 -></i - wi 4) \h-) = - i^ 4) ifa-><fa - i^i/t-) 

l ' 2 ~\(2hp 3 ) (2/2*,) + (2/ 2jP3 ) (2/i*;) -(2/1/2) (2p 3 k l )]. (48) 



(/2-|^ 3 |/l- 



For the first term one uses in turn again eq. (42), while the last term in the square bracket reduces 
the rank by one through eq. (15) and eq. (29). This allows to reduce any rank r > 2 integral to 
scalar or rank 1 integrals. It remains to treat rank 1 integrals. For rank 1 integrals we may use 

< ' 2 - |</ ' (4, "'- > = ik-^H-M ^f^)' m 

where the subscript "+" indicates that a projection operator (1 + Ys)/2 has been inserted into the 
trace. Since the piece proportional to the totally antisymmetric tensor vanishes after integration, 
we may replace 



Tr + (V! 4 W 3 f 1) - ^Tr (VJ'Wsf 1) , Tr + ty%M 2 ) - ^Tr (V z (4) f i// 3 f 2) ■ 



(50) 



Therefore 



/(4)|- 1 



(h - WT%-) = ofi \^T~\ [(^iP3)(2/ 2 */) + (2/2/73) (2/i*/) - (2/1/2) (2p 3 */)] 

2(/2- I/73K1-) 



+ terms, which vanish after integration, 
1 

2(h-W 3 \l 2 -) 

+ terms, which vanish after integration. (51) 



(h ~ |#! 4) |/i-> = ^ 1, , [(2/ip 3 ) (2/2*0 + (2/ 2j p 3 )(2/i*/) - (2/1/2) (2p 3 *0] 



This allows to reduce rank 1 four-point integrals. 



3.5 The tensor five-point function and beyond 

Here we discuss the tensor reduction of n-point functions with n > 5. For rank r > 2 we follow 
the same steps in eq. (47) - eq. (48) as for the four-point function. The only difference occurs in 
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the treatment of rank one integrals. We note that for the n-point functions with n>5we have one 
further additional independent momentum, which will be labelled p$. For the rank one integrals 
we have [29] 

(h - \Hf ] \h-) = ~\ [(2/1/74) (2l 2 ki) + (2l 2 p 4 ) (2hk) - (2hl 2 ) (2p 4 k l )} (h - \rf3\h-) 
+^ [(2/1/73) (2/z*,) + (2/ 2jP3 ) (2Zi*,) - (2hl 2 ) (2p 3 k l )} (h - W4\h~), 

(h ~ W^lh-) = I [{2h P A) (2hki) + (2/ 2jP4 ) (2hk) - (2/ife) (2 P 4k)] (h - W3\h~) 

~ [(2/i/7 3 ) (2hki) + (2/2/73) (2/i*,) - (2/1/2) (2/7 3 */)] (/2 - I//4I/1-), (52) 

where 

5 = (/i-|// 4 |/2-)(/2-|// 3 |/l-)-(/l-|//3|/2-)(/2-|//4|/l-). (53) 

8 is proportional to the square root of the Gram determinant of the four-momenta l\ , I2, P3 and p$. 
Numerical instabilities in the limit 8^0 can be treated with the methods discussed in ref. [29]. 

4 Higher-dimensional integrals 

In this section we discuss the evaluation of scalar integrals of the form 

iDu f-k 2 



r0,s _ zj E 2e f ±_k \ t^h m , 

" ^ J * 2 (*-pi) 2 ...(*-pi-...p„-l) 2 ' P J 

with s > 0. Scalar integrals with s = are treated in section 5. In a space of D = 2m — 2e 
dimensions (with m being an integer), we decompose k 2 D ^ as follows: 

rfp) = ^(2m)+^(-2£) (55) 

If a power of (— * 2 _ 2e )) appears in the numerator we have [58] 

rd 2m -^k. , 2 , Sf/t u , 2 N r(j-e) r d 2m+2s - 2£ k fnJI , 2 x , c ~ 

y ^r ( " fc (- 2e ) ) = TFiT J * m+s - e i /( ^)'^- 2£ ) ) - (56) 

The effect of a factor of (— * 2 _ 2e p s m the numerator is to shift the dimension by 25. Note 
that r(s — e) /r(— e) brings an explicit factor of £, therefore we have to take higher-dimensional 
integrals into account only if they are divergent. A scalar n-point integral with unit powers of the 
propagators is finite, if [27] 

2 < I < n. (57) 
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Here 2 < D/2 is the condition to be infrared finite and D/2 <nis the condition to be UV-finite. 
Therefore, higher dimensional integrals are always infrared finite and we only have to calculate 
the UV-pole of the higher dimensional integrals. This can easily be done. For m>nwe find 



d 2m ~ 2e k 1 



in m ~ e k\k\...kl 

1 (-l) m 
£ (m — n)\ 







j d n ad\ 









(58) 



where 



7 = -Y, a i a j(Pi+i + -+Pj) 2 ■ (59) 

i<j 

Note that the integral over the Feynman parameters is a polynomial in the Feynman parameters 
and can be done according to the formula 



/ 




/flSli-Vfl.lfl] 1 - 1 ...^- 1 = r (vi)-r(v„) 
^ ' ' 1 r(vi + ...+v M ) 



In practice there are additional simplifications: When calculating one-loop amplitudes, we are 
free to choose an appropriate gauge. Using the Feynman gauge, we can ensure that the rank r 
of a loop integral is always less or equal the number of external legs n. In addition, there are 
obviously no powers of ^_ 2e ) m the original loop integral, i.e. we have 

r<n and 5 = 0. (61) 

The algorithm for the tensor reduction in section 3 respects the inequality 

r + 2s<n. (62) 

Therefore the only non-zero higher-dimensional integrals which occur in the Feynman gauge 
result from the two-point function with a single power of ^_ 2 e) m the numerator (n = 2 and 

5=1), the three-point function with a single power of &(_ 2e ) m the numerator (n = 3 and 5=1) 
and the four-point function with two powers of ^_ 2 e) m the numerator (n = 4 and 5 = 2). The case 
of the two-point function has already been discussed explicitly in section 3.2. For the remaining 
two cases one finds: 

, ,iD h ( -k 



£ V £ I u K ^"^^ = - + 0(8), 
J z7t§ kfklkj 2 

d D k (-*( 
z7t"2 k^k^k^k^ 
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5 Reduction of higher point scalar integrals 



In this section we discuss the reduction of scalar integrals of the form 



d D k 1 
-q k Hk- Pl )i...(k-pi-...Pn-i) 2 ' 



I n = e^^j d -4— (W) 



with n > 5 to a basic set of scalar two-, three- and four-point functions. It is a long known fact, 
that higher point scalar integrals can be expressed in terms of this basic set [59, 60], however the 
practical implementation within dimensional regularization was only worked out recently [25, 
37^4-0]. We distinguish three different cases: Scalar pentagons (i.e. scalar five-point functions), 
scalar hexagons (scalar six-point functions) and scalar integrals with more than six propagators. 



5.1 Reduction of pentagons 

A five-point function in D = 4 — 2e dimensions can be expressed as a sum of four-point functions, 
where one propagator is removed, plus a five-point function in 6 — 2e dimensions [37]. Since the 
(6 — 2e) -dimensional pentagon is finite and comes with an extra factor of £ in front, it does not 
contribute at (e°) . In detail we have 

/ 5 = -2eB/ 5 6 - 2e + £^/f = £ Vf + 0(e), (65) 
i=\ i=i 

where /f _2e denotes the (6 — 2e) -dimensional pentagon and 1^ denotes the four-point function, 
which is obtained from the pentagon by removing propagator i. The coefficients B and bi are 
obtained from the kinematical matrix Sij as follows: 

b i = %( S ~ 1 )i j > B = L b - ( 66 ) 

j i 

5.2 Reduction of hexagons 

The six-point function can be expressed as a sum of five-point functions [38] 

6 n 

h = i>/r. (67) 

The coefficients bi are again related to the kinematical matrix S^: 

bi = L(S-%- (68) 

j 
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5.3 Reduction of scalar integrals with more than six propagators 

For the seven-point function and beyond we can again express the n-point function as a sum over 
(n — 1) -point functions [40]: 

In = t^n-V ( 69 ) 

i=l 

In contrast to eq. (67), the decomposition in eq. (69) is no longer unique. A possible set of 
coefficients r ; can be obtained from the singular value decomposition of the (n — 1) x (n — 1) 
Gram matrix 

4 

Gij = Y, U ^k{V T ) (70) 

k=\ 



as follows [27] 



with 



r,- = — , 1 < i < n — 1, 

wy — 

r n = (71) 

7=1 



^5 = ^GjjVjs. (72) 



Note that the kernel of G ;J - is spanned by the vectors V^, V^, Vu n _iy 



6 Numerical implementation 

We have implemented the algorithms described so far into a numerical computer program. The 
program is able to calculate the coefficients C_2, C_i and Co of the Laurent expansion of one- 
loop n-point integrals of rank r and s powers of ^_ 2 e) m the numerator: 



17 = ^ + ^-+C + o(£). (73) 



As our algorithms are valid for any number of n external particles, the actual limitation on n will 
result from the available computer power. 

We have performed several checks on our computer code. The value of a tensor integral 



^( a ,- fel | fel -)...(.-ivi^-)/g ^_ ; , i)2 5j; _,„ p „_ i)2 
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is clearly unchanged if we permute the tensor structure 

(a\ ~ |Y W I bi-) ... (a r - \j Mr \b r -) -> (o o(1) - |y w | b a[l) -) ... (a a[r) - \% r \b o{r) -) . (74) 

Since our algorithm reduces the rank step by step, this actually provides a non-trivial check. 
Secondly, for specific choices of the tensor structure, like 

{ Pi -\Hj\pi-) = 2 Pi kj, (75) 

the numerator reduces immediately to simpler integrals. This will lead to relations among differ- 
ent integrals, which can be checked numerically. 

Finally, we have written three independent codes (in two different programming languages: 
Fortran and C++), which all agree with each other. 

For future reference we give a few numerical results. We start by specifying a set of twelve 
light- like momenta p,, with i = 1, 12. This serves as the input data for the scalar 12-point 
function, where all external particles are light-like. By combining four-vectors we can obtain the 
external kinematics of lower point functions. We choose the set 

{ph-,Pj-uPj + - + Pn} (76) 

for the j-point function. For 3 < j < 12 this corresponds to j — 1 light- like external legs and one 
massive leg. Note that for j = 2 we have a two-point function with light- like external momenta, 
which vanishes. The random values for our set of momenta (in units of GeV) are: 

pi = (5.897009121257959,-1.971772490149703,-4.63646682189329,-3.064311543033953), 

p 2 = (9.78288114803946,-3.495678805323657,-7.42828599660035,-5.320297021726135), 

P3 = (3.751716626791747,0.3633444560526895,-2.74701214525531,2.529285023049251), 

p 4 = (14.8572007649265,-9.282840702083684,9.182091233148681,7.08903968497886), 

p 5 = (4.056006277332882, -1.236594041315223,0.8781947326421281,3.761754392608195), 

p 6 = (2.023022829577847,0.3217130479853592,0.6516721562887716, 1.887973909901054), 

Pl = (23.51469894530697,20.57030957903025,3.67304549050126,-10.78481187299925), 

P& = (6.161822860155142,0.9716060205020823,2.082735149413637,5.717189606638176), 

p 9 = (10.67981737238498, -2.237405231711613,-2.487945529176884,-10.14212226215274), 

pio = (9.275824054226526,-4.002681832986503,0.83197173093136,8.326300082736525), 

pn = (-45,0,0,45), 

pn = (-45,0,0,-45). (77) 
This set satisfies momentum conservation 

12 

J> y - = 0. (78) 

7=1 

We give values for tensor integrals up to rank 2. Since for higher rank integrals no new reduction 
algorithms are used, this is sufficient for demonstration purposes. For j-point integrals with 
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j < 10 the momenta p\\ and p\i have no special relation to the external kinematics (only the 
sum pj + ... +pn corresponds to an external leg). Therefore the sandwich 

(Pi2-Wi\pn-) (79) 
is an example of a generic rank 1 integral. Similar, we use for j < 8 the tensor structure 

(pi2-Wi\pn-)(pw-Wi\P9-)- (80) 

The numerical values of the bra- and ket-spinors depends on a choice for the phases of the 
spinors. Our conventions are listed in the appendix. In addition we use, when evaluating spinors, 
a rotation (jc, y,z) — > (z,x,y) for the spatial coordinates of a four- vector, such that the line, where 
spinors are not defined, lies along the negative y-axis. This avoids problems with incoming 
particles, which are often taken to be on the z-axis. For cross-checks we also quote the numerical 
values of the spinors in our convention: 

( Pn -\ = (-6.708203932499369,-6.708203932499369), 

( Pw -\ = (3.179275984427568,2.618929631626706 + 1.258991623436299/), 

\ pn -) = (6.708203932499369,-6.708203932499369), 

\ P9 -} = (2.862144623041976,-3.543539407653475-0.7817233321122782/). (81) 

The results of the loop integrals will depend also on the renormalization scale (jl. We set 

/j = 135 GeV. (82) 

This specifies all input parameters. The results for the coefficients C_2, C_i and Co of the Laurent 
expansion are shown for the scalar integrals in table 1 . The corresponding numbers for the rank 
1 integrals can be found in table 2, while table 3 shows the results for the rank 2 integrals. Our 
independent programs agree within 10~ 7 . Table 4 shows the CPU time in seconds for a tensor 
integral with n external legs and rank r for r < n < 10 on a standard PC equipped with a Pentium 
IV running at 2 GHz. The recursive algorithm is efficiently implemented with the help of look-up 
tables. The required memory for the look-up tables is negligible, i.e. of the order of 10 MB for 
the case n = r = 10. 



7 Conclusions 

In this paper we discussed an algorithm for the automated computation of one-loop integrals, 
which occur in a massless quantum field theory. This is relevant for high-energy experiments, 
where the masses of the quarks (with the exception of the top quark) can usually be neglected. 
We reported on the implementation of this algorithm into a numerical program. It is worth 
to point out, that there are a priori no restrictions on the number of external legs of the loop 
integrals. Therefore the actual restriction is only given by the available computer resources. We 
gave examples for the evaluation of loop integrals with up to twelve external legs. In future 
work we intend to integrate this program into a package for the automatic calculation of jet cross 
sections. 
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n 



3 

4 

5 

6 

7 

8 

9 

10 

11 

12 



C 7 



9A321 ■ 10° 
3.0405-10 1 
1.3863 • 10 3 
6.7736 -10- 6 
3.1190- 10 8 

1.68789- 10 11 

9.56014 

4.41232 

1.51680 
-8.17311 



io- 15 
io- 18 
io- 21 
io- 25 



C_i 

(1.1371 +0.2963/) -IO 2 
(3.8038 + 0.9552/) -10° 
(1.8815 + 0.4355/) • 10~ 2 
(9.8029 + 2.1280/) -IO" 5 
(4.6404 + 0.9799/) -IO" 7 
(2.4470 + 0.5303/) -10 10 
(1.2782 + 0.3003/) • 10" 13 
(5.3676 + 1.3862/) -IO" 17 
(1.7511+0.4765/) - 10 20 
(-9.3478- 2.5677/) -IO" 24 



Co 



~ (6.3106 + 3.5723/) • 10 2 
(2.2164+ 1.1950/) -IO 1 
(1.2211 +0.5911/) -IO" 1 
(6.8667 + 3.0797/) -IO" 4 
(3.3624+ 1.4578/) • IO" 6 
(1.7380 + 0.7688/) • 10~ 9 
(8.4151 +4.0157/) -IO" 13 
(3.1673 + 1.6863/) - 10 16 
(9.6086 + 5.5011/) • IO" 20 
(-5.0527- 2.9367/) -IO" 23 



Table 1: Results for the scalar n-point functions with 3 < n < 12. The Q denote the coefficients 
of the Laurent series. 



3 
4 
5 
6 
7 
8 
9 
10 



C-2 C l Co 

(2.3701 + 1.2937/) • 10 3 (2.9247 + 2.5629/) • 10 4 

(-1.0164 -0.4783/) • 10 2 (-1.1710 -0.9390/) • 10 3 (-6.1191 -7.9185/) • 10 3 

(-7.6639 - 3.6047/) • IO" 1 (-9.7028 - 7.5033/) • 10° (-5.7693 - 6.8698/) • 10 1 

(-4.3386 - 2.0356/) • IO" 3 (-5.7717 - 4.3739/) • 10~ 2 (-3.6478 - 4.1735/) • IO" 1 

(-2.0606 -0.9984/) • IO" 5 (-2.7898 - 2.1356/) • IO" 4 (-1.8100 - 2.0589/) • IO" 3 

(-1.0155 -0.6865/) • 10~ 8 (—1.3143 - 1.2571/) • 10~ 7 (-0.8305 - 1.1134/) • 10~ 6 

(-4.0783 -6.1625/) • IO" 12 (-4.1976 - 8.6804/) • IO -11 (-2.1879 -6.2292/) • 10~ 10 

(-1.0266 - 3.9800/) ■ 10~ 15 (-0.2878 - 4.7984/) ■ 10~ 14 (0.2754 - 2.8456/) ■ 10~ 13 



Table 2: Results for the n-point functions of rank 1 . The Q denote the coefficients of the Laurent 
series. 



C-2 C-i Co 

(3.1317 + 4.3445/) • IO 3 (3.3251 +7.4740/) • 10 b 

(-0.7042- 1.0292/) • 10 4 (-0.6161 - 1.5979/) • 10 5 (-0.1596 - 1.1662/) • 10 6 

(-5.3188 -7.7592/) • 10 1 (-0.5233 - 1.2881/) • 10 3 (-0.1818 - 1.0225/) • 10 4 

(-3.0368 - 4.3882/) • IO" 1 (-3.1500 - 7.5262/) • 10° (-1.2159 - 6.2003/) • 10 1 

(-1.4611 - 2.0679/) • 10~ 3 (-1.5537 - 3.6078/) • 10 2 (-0.6331 -3.0290/) • IO" 1 

(-6.5351 -9.7049/) • IO" 7 (-0.7106 - 1.7046/) • 10~ 5 (-0.2959 - 1.4465/) • IO" 4 



Table 3: Results for the n-point functions of rank 2. The Q denote the coefficients of the Laurent 
series. 
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x 

n 







1 






2 




3 




4 


5 


6 


7 


8 


9 


10 


2 


1 • 10 


-6 


5 


10 


-6 


2- 10" 


-5 




















3 


1 • 10 


-6 


2 


10 


-5 


2- 10" 


-4 


1 • 10 


-3 
















4 


2- 10" 


-6 


5 


10 


-5 


4-10" 

"T 1-\J 


-4 


2- 10" 


-3 


6 • 10 3 














5 


3 • 10" 


-5 


1 


10 


-4 


6-10 


-4 


3 • 10" 


-3 


9- 10~ 3 


0.03 












6 


2-10" 


-4 


3 


10 


-4 


9-10 


-4 


4-10" 


-3 


0.02 


0.04 


0.1 










7 


7-10" 


-4 


7 


10 


-4 


1 • 10 


-3 


5-10" 


-3 


0.02 


0.06 


0.1 


0.4 








8 


3 • 10 


-3 


3 


10 


-3 


4-10" 


-3 


8 • 10 


-3 


0.02 


0.07 


0.2 


0.6 


1.8 






9 


0.01 




0.01 




0.01 




0.02 




0.03 


0.08 


0.3 


0.9 


2.6 


7 




10 


0.05 




0.05 




0.05 




0.06 




0.06 


0.2 


0.4 


1.1 


3.5 


8 
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Table 4: CPU time in seconds for a tensor integral with n external legs and rank r for r < n on a 
standard PC ( Pentium IV with 2 GHz). 



A Complex four- vectors and spinors 

In this appendix we list our conventions for spinors. In particular we comment on complex four- 
vectors and their associated spinors. Although the external momenta of a loop integral are real 
quantities, the decomposition of two massive vectors into linear combinations of null-vectors, as 
in eq. (24), may introduce complex four- vectors. For the metric we use 

8lN = diag(+l, -1,-1,-1). (83) 

A null-vector satisfies 

(po) 2 -(pi) 2 -(P2) 2 -(P3) 2 = 0. (84) 
This relation holds also for complex p^. Light-cone coordinates are as follows: 

P+=Po + P3, P-=Po~P3, P±=Pi + ipi, P±* = pi-ipi- (85) 

Note that p±* does not involve a complex conjugation of p\ or p2- We use the Weyl representa- 
tion for the Dirac matrices 

v=(£7), *=*VtV=(; »). 

where the 4-dimensional G^-matrices are 

o^=(l,-6), a^ s = (l,6), (87) 
and Pauli matrices a = (a x , a y , a z ) are as usual 
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Four-component Dirac spinors are constructed out of two Weyl spinors: 

(£!)-(?)-(:$)■ 



where 



u±{p) = 1 -(\±y 5 )u(p). (90) 

Bra-spinors are given by 

u(p) = ({p-\AP+\)= (p A ,Pb) = (u-(p),u+(p)), (91) 

where 

u±(p) = U(p) l -(\Ty5). (92) 

Eq. (89) and eq. (91) show three different notations for Weyl spinors. We are using mainly 
the bra-ket notation. In terms of the light-cone components of a null-vector, the corresponding 
spinors can be chosen as 

b+> = ' ( 1 , , p _; - ( p* 



^/\p~A \ p+ J' V\pA V P± J ' 

(P+\ = -7j=f (-P±,P+) J = ^^(P+,PJ-*) 5 (93) 

vb+l vIp+I 

where the phase § is given by 

P+ = \p+\S*- (94) 

The spinor products are then given by 

(pq) = (p-\ q +)= 1 -A p±*q+~P+q±*), 

v\p+\ k+l 

M = («+|p-> = - 7 =i=e-'*'e-^(pj.« + -p + «j.). (95) 

B The basic scalar integrals 

In this appendix we list the basic scalar integrals, which are the scalar two-point, the scalar three- 
point and the scalar four-point functions in D = 4 — 2e dimensions. Since we restrict ourselves 
to massless quantum field theories, all internal propagators are massless and we only have to 
distinguish the masses of the external momenta. All scalar integrals have been known for a long 
time in the literature. Classical papers on scalar integrals are [61, 62]. Scalar integrals within 
dimensional regularization are treated in [37,63]. Useful information on the three-mass triangle 
can be found in [64-66]. The scalar boxes have been recalculated in [67, 68]. 
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B.l The two-point function 

The scalar two-point function is given by 



HpIi?) = i + 2-ln 



1 



(96) 



B.2 Three-point functions 

For the three-point functions we have three different cases: One external mass, two external 
masses and three external masses. The one-mass scalar triangle with p\ ^ 0, p\ = p\ = is 
given by 



z L p\ epf V P 



p 2 A + ± ln 2(-p* 



2pj 



1 



2p 2 



+ 0(e). (97) 



The two-mass scalar triangle with p\ ^ 0, p\ ^ and p\ = is given by 

1 1 



j-2m/„2 2 2\ 



Hp 2 -Pi) 
l 



2\ / 2\ 1 

p A , / -Pi 



M+ln 



-In 



+ 0(8). 



(98) 



The three-mass scalar triangle with p\ ^ 0, p\ ^ and ^ 0: This integral is finite and we have 



^{pIpIpW) = -Jd 3 a- 



8(l-ai-a2-a 3 ) | p(g) 
aia 2 pj -a 2 a 3 ^-a 3 ai^ 



(99) 



With the notation 

Si /?|, h = pj-pi-p 2 , h=p 2 -p 2 -ph 

A 3 = W) 2 + {plf + {Pi) 2 - 2p\pl - 2p 2 p\ - 2p 2 3 pl 
the three-mass triangle /| ra is expressed in the region p\ ,p\,p\ < and A 3 < by 
2 



j3m 



X 



CI2 I 2 arctan 



5i 



+ C1 2 ^2 arctan (^^)) + C1 2 (2 



(100) 



arctan 



5 3 



+o(e). 



(101) 



The Clausen function Cl2(x) is defined in eq. (119). In the region p\ ,p\,p\ < and A 3 > as 
well as in the region p\,p 2 <0,p\>0 (for which A 3 is always positive) the integral /| m is given 
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by 



j3m 



Re 



2(Li 2 (-px) + Li 2 (-p V )) + ln(px)ln(Py) + h ( ^ ^ + 3 1 



■In 



A 3 1(5 1 - v ^)(5 3 -v / a;) 



+ o(e), 



where 



2p 2 



2 2 
Pi ^2 n 

pf /7| 8 3 + VA 3 



The step function 0(x) is defined as 0(x) = 1 for * > and Q(x) =0 otherwise. 



(102) 



(103) 



B.3 Four-point functions 

For the four-point function we use the invariants 



2 2 

s= (P1+P2) , t= (P2+P3) 



together with the external masses m? = pf. 
The zero-mass box (m 2 = = m 2 = m\ = 0) 

4 2 



/ 4 °-(s,^ 2 ) = 



£ 2 5? est 



>"(?) +ln (7 



1 

H — 
st 



—5 



The one-mass box (m 2 = = m 2 = 0): 



*4 {sAmW) = ^--- 



e z st est 



ln( )+ln( -4 ) -In 



1 

+ - 

5? 



(104) 



+ o(e). (105) 



+ ln 2 



-In 



2 / 1 2 
/ VP 



2 / m 4\_ | n 2 ' ' v 



-3^ + 0(8). 
The easy two-mass box (m 2 = m 2 = 0): 

2 



(106) 



/ 2rae (M,m 2 ,m 2 ^ 2 ) 



1 



5? — m^rn 2 ^ 



-2Li 2 1 



In 2 



In 2 



2" i 4; 



In 



-m; 



-In 



In 2 



— m% 



In 2 



(-*) ; v (-0 

2Ti (\ ( " m 4^ + 2Ti (\ (- m 2)(" m 4 



(-m 2 4 ) 



( - m2 ^-2Li 2 fl-t^)-2Li 2 (l-^i 



+ 0(8). 



(107) 
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The hard two-mass box (m\ = m\ = Q): 

1 1 



Il mh (s,t,mlmW) 



z L st est 



ln(^-]+21n(-| 



1 

+ - 

st 



-hr( — )+hr(^2 



In 



In 



-m; 



A< z 



A< z 



1 "7 



ln( -4 )ln 



-2Li 2 ( 1 



-m. 



{-ml) 
i-t) 



The three-mass box (m^ = 0): 
/ 4 3m (s,t,n%,n%,ml,i?) =- 



(-^) 



2Li 2 1 



(-0 



2^2 



+ o(e). 



(108) 



£ [st — rn^rnfj 



In ^ +ln ^ 



st — m^rn^ L 
!„(-?] In 



—s 



-\n~l- - )+-ln-( ^2 



-In 2 



In 



In 



-m; 



In | -2 I in 



A< z 



+ ln 



In 



-ln z 



— m; 



-In 



—s 



In- ( — 



A< z 



In 



^ z 



(-*) 



— m\ 



( (-mi)\ ( (-mi)\ ( (-ml) (-ml) 



+ 0(e). (109) 



(-0 J V (-*) (-0 

The four-mass box: 

(5,?,m2,m3,m4,/j 2 ) = I^ m {st,m\m^ : m^m^,^') +K{s,t,m\,m^,m\,m^), (110) 

where 



2jii 



: 3 



#(si,*l,S2,'2,53,'3) = -^-£9(-^)0Hi) 



(=1 



X 



ln £ - - X) ( 1 + 10) - In £ Sjtj - ( Si ti + l)(l+ iO) 



(HI) 



and 



X = \J (s\t\) 2 + (S2t2) 2 + (S3t3) 2 - 2s\t\S2tl ~ 2s2^3?3 ~ 253^3*1^1 • (H2) 
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C Analytic continuation 

In one-loop integrals the functions 



ln( — f 



Li 2 1 



(zf) 
(-0 



(113) 



and generalizations thereof occur. The analytic continuation is defined by giving all quantities a 
small imaginary part, e.g. 



s — > s + iO. 

Explicitly, the imaginary parts of the logarithm and the dilogarithm are given by 



(114) 



In 



Li 2 1 



—s 

-t 
s) 



= In 



(-0 

This generalizes as follows 

ln A-*i)(-*2) 



)-to[&(s)-B(t)], 
ReL l2 (l-9-,-e(-9ln(l-9lmln(^ 



(115) 



In 



SlS 2 



ht 2 



(-h) i-ti) 

V (-h)(-t2)J 2 \ ht 2 ) V (-h)(-h)J V 



-iK[8(si) + e(s 2 )-e(ti)-8(t2)], 

S\S 2 



(-*) tf2) 

ho (-fe)y 



where 



In 1- 



(~*l) (~^2) 
Hi) ("fe) 



) 


1 S1S 2 







- \in [0(50 + Q(s 2 ) - 6(fi) - 6(fe)] 6 ( ^ - 1 



D Numerical evaluation of special functions 

The real part of the dilogarithm Li 2 (x) is numerically evaluated as follows: Using the relations 



Li 2 (x) = -Li 2 (l — x) + — -ln(x)ln(l —x), 

b 

= - Li 2(£)-y4 (ln( ~* ))2 ' 



(116) 



the argument is shifted into the range —1 <x< 1/2. Then 

Bi 



Li 2 (x) = £ 



1 -z l+l 



2 - + ^(2n+l)!' 



(117) 
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with z = — ln(l — x) and the Bi are the Bernoulli numbers. The Bernoulli numbers 5, are defined 
through the generating function 



t t n 



It is also convenient to use the Clausen function C\2{x) as an auxiliary function. The Clausen 
function is given in terms of dilogarithms by 

Cl 2 (x) = i [02(^-02 (*-*)]. (119) 
Alternative definitions for the Clausen function are 

X 

Cl 2 (x) = £ = - J dtln (|2sin (0 |) . (120) 

n=l Q 

The Clausen function is evaluated numerically as follows: Using the symmetry 

Cl 2 (-x) = -CI2W, (121) 

the periodicity 

Cl 2 (x + 2rm) = Cl 2 (x), (122) 

and the duplication formula 

Cl 2 (2x) = 2Cl 2 (x)-2Cl 2 (7t-x) (123) 
the argument may be shifted into the range < x < 2n/3. Then 

C1 2W = -,,„ W+ , + |y^» + .. (124) 
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